Tuning impurity states in bilayer graphene. 
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We study the impurity states in bilayer graphene in the unitary hmit using Green's function 
method. Unlike in single layer graphene, the presence of impurities at two non-equivalent sites in 
bilayer graphene produce different impurity states which is understood as the change in the band 
structure due to interlayer hopping of electrons. The impurity states can also be tuned by changing 
the band structure of bilayer grahene through external electric field bias. 

PACS numbers: 



I. INTRODUCTION 

In recent years, the fabrication of few layers of 
graphene systemsii^i^ has attracted a lot of attention to 
study the electronic properties of these systems. The 
electrons in single graphene show some unconventional 
electronic properties such as the half integer quantum 
Hall effect, and Klein paradox. The electrons in bilayer 
systems also show some interesting properties as seen in 
quantum Hall effect. Many properties of single and bi- 
layer graphenes also differ because of the different crystal 
structure. In this communication we present a systematic 
study of impurity effects in bilayer graphene and contrast 
it with those of the single layer counterpart. 

The fundamental difference between single and bilayer 
graphene originates from their crystal structure. Sin- 
gle layer graphene is an atomically thin two-dimensional 
hexagonal packing of sp^ bonded carbon atoms. It is the 
building block of (multi-layer) graphene. One unit cell 
of single layer graphene has two non-equivalent lattice 
sites (A, B). As s result, the electron wave function is 
spinor-like, where the sublattice index plays the role of 
pseudo-spin. The tight binding calculation^ shows that 
the electrons in single layer graphene disperse linearly, 
i.e. {Ek = ±VFk), where vp = ^ = 5.8 eVA is the 
Fermi velocity, t = S.OeV is the nearest neighbor hop- 
ping energy in the plane, and a is the lattice constant, 
and hence are called the massless Dirac fermions. The 
bilayer graphene, as shown schematically in Fig. [l] can 
be thought of as a stacking of two identical single layer 
graphenes in the third dimension. In one of the common 
ways of layer stacking, known as Bernal stacking, only 
one of the non-equivalent lattice sites (site-A) stay on 
top of each other, another site (B) lies in the middle of 
the hexagon of the other layer The electron can hop be- 
tween the layers along the bonding of these two A-sties 
with a hopping energy {tj_), which is about ten times 
smaller than the hopping energy along the plane. This 
interlayer hopping hybridizes the pz orbital of the carbon 
atom at site-A resulting in different dispersion relation 

of the electrons, {Ek ^ ± \J + v^k"^) Two 
of the branches of the electronic band touch each other 
at the Fermi energy, whereas the other two branches be- 
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FIG. 1: (Color online) Lattice structure of bilayer graphene 
under consideration for a tight binding calculation. The green 
(dashed) line forms the top layer and the red (dotted) line 
forms the bottom layer. 



come gapped with energy gap equal to t^. The disper- 
sion relation of the electrons corresponding to the gapless 
branches can be expressed in parabolic form at low mo- 
menta. The electron energy in this case have parabolic 
dispersion at low momenta. In addition to this differ- 
ence of the band structure between single and bilayer 
graphenes, the bilayer system gives a freedom of tailor- 
ing the band structure by applying an external electric 
field bias (V) on the two layers. The dispersion relation 
in the presence of the external field bias becomes^ 



(1) 

which shows that the two valence (conduction) bands 
shift below (above) the Fermi energy by V/2. 

To understand the electronic property of these systems, 
it is important to study the revealing impurity effects. 
The impurity states in single layer graphene has been 
studied by Wehling et al}^ and by Bena.^^ It has been 
shown that if an impurity is introduced on a lattice site, 
a virtually bound impurity resonance state can be in- 
duced, which can be seen as an enhancement of the local 
density of states (LDOS) at the neighboring site. This 
effect is symmetric with respect to the sublattice the im- 
purity is on, i.e., an impurity at site B shows exactly the 
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same effect as it is at site A. In bilayer graphene, due 
to the bonding between the sites A of the two layers, 
the site A and B in each layer are no longer equivalent. 
The interlayer hopping distinctly differentiates the two 
sites. So we expect to see different impurity states with 
respect to the location of the impurity. In addition, it is 
natural to expect that the external bias which changes 
the band structure also modifies the properties of the im- 
purity states. This motivates us to study the impurity 
states in bilayer graphene as a function of both t± and V. 
The impurity effect in single layer graphene is contained 
in our discussion as a special case of bilayer graphene 
where t± = 0. 

The outline of the paper is as follows: In Sec. II, we 
present a Green's function method to study the impu- 
rity effects. For the clean system, the results from the 
Green's function method is benchmarked with the exact 
diagonalization. In Sec. Ill, we present a systematic set 
of results and discussions. Concluding remarks are given 
in Sec. IV. 



using Green's function. We need to find the form of 
the Hamiltonian to define the free particle Green's func- 
tion. The Green's function in the presence of the impu- 
rity potential is obtained by using the Dyson's equation. 
We derive a tight binding Hamiltonian using the lattice 
structure of the bilayer graphene as shown in Fig. [T] 
The basis vectors are chosen to be, ai = (^, ~^° ), 
and fli = i^,^^)- The correspond basis vectors of 
the reciprocal space is given by bi = ^^^d 

62 = (|^'"7f^)- Sites Ai and A2 are connected along 
the z-direction, the electrons have t± hopping energy in 
this direction. We have assumed that the layer which 
has lattice sites Ai, Bi is biased with V/2 and the layer 
which has lattice sites ^2, B2 with —V/2 so that the po- 
tential difference between the two layers is V . For the 
convenience of discussion, these two layers are called the 
top and bottom layers, respectively. 



II. THEORETICAL METHOD 



In this section we describe the theoretical method 
used to calculated the LDOS in the bilayer graphene 



The tight binding equations for this lattice structures 
in the clean case can be written as. 



J 



— (/)Ai -i[2exp( — ^)cos(^^-^ — ) -I- exp(ifc^a)](^_Bi - t_L<^A2 =-E0Ai, 



(2a) 



V 



- t[2 exp(— — ) cos(^^- — ) + exp(-ifc;j;a)]0Ai + tt'^Si = E(f>B^ 



(2b) 



■ t±(pAi + -7r^A2 - R2exp(— — )cos(-^^— — ) 



-I- exp(-ifc^a)](/)B2 = E(i>A2, 



(2c) 



,r„ , -ikxa kya\/3 
t[2exp(^— )cos( 



2 / + eyip{ikxa)](t)A2 + { — )<pB2 = E(j)B.2 ■ 



(2d) 



These equations can be written in a matrix form. 



Here, 
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(3) 



with t = -t[2exp(^^)cos(^^^) -h exp(iA::ra)] and t* 
is the complex conjugate of t. 

The Green's function corresponding to this Hamilto- 
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nian can be expressed as, 



(a) U=0, tj^=0.0, V=0.0, DOS at all sites (b) U=0, t^-O.O, V=0.5, DOS in two layers 



G(°)(k,L.) = M-i?k]- 



(5) 



where 1 represents the unit matrix. As long as aU in- 
troduced impurities (up to 4) are located within a single 
cell, we can express the local Green's function exactly 
through the T-matrix method, which leads to 



(c) = Gf (c) + Glo"^ (c.)T»G°^. (c.) 



(0), 



(6) 



Here the T-matrix f (w) = (7(1 - G^°'> {uj)U)-^ , the lo- 
cal Green's function gI^\u!) = jjJ2w^^^^0^^^) where 
the summation is over the first Brillouin zone, and U 
is the matrix representation of the impurity potential. 
The local density of states at different sites is given 
by, Nai = ^Gii(w -I- ij), Nbi = ^G22(t^ + h), 
Na2 = ^G:i:i{uj + 17), Nb2 = ^Gii{uj + ij), where 
7 is the lifetime broadening. For numerical calculation 
all relevant energies are measured in terms of the in- 
tralayer hopping energy of the electron, t. For simplicity 
we use the impurity potential close to the unitary limit, 
U = 100. The number of 1024 x 1024 k points are used 
in the Brillouin zone. The intrinsic lifetime broadening 
of 7 = 0.005 is taken. When we present results, unless 
otherwise stated, the panel a to d represent the LDOS 
at sites Ai, Bi, A2, and B2, respectively. In each figure, 
the DOS is plotted (along the vertical axis) as a function 
of energy E (along the horizontal axis) . 

We can also calculate DOS in the absence of impu- 
rities using eigenvalues and eigenvectors of Hamiltonian 
matrix, Eq. A little algebra yields the eigenvalues 



(7) 



where £„ represents the four eigenvalues. Note that the 
linearization of Eq. ([7]) near the corner of the Brillouin 
zone reduces to Eq. ([1]). 

The corresponding eigenvectors are given by a general 
equation. 





/ 



I? 



\ 2 ^ ^fi 



V 



r+£n 
t 
1 



(8) 

where substituting four different values of e„ gives the 
four sets of eigenvectors. We normalize thus obtained 
4x4 matrix of the eigenvectors and calculate the density 
of states on four nonequivalent sites using the standard 
equation. 



Ntt 
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(ii;-£„)2+72 



(9) 



where a — A, B, P = 1,2, (pap^e^lli) is the normalized 
eigen- vectors. We used this DOS to check the result ob- 
tained by Green's functions method in the absence of 
impurity. 
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FIG. 2: (Color online) The density of states as a function of 
energy is shown, a) For t_L = the DOS at all sites is equal, 
and vanishes linearly at the Fermi energy, b) For t± — 0, 
and V 0, the DOS curve for sites of the top (bottom) layer 
shifts above (below) the Fermi energy by V/2 creating finite 
DOS at the Fermi level, c) For tx / DOS at A sites of both 
layers is equal and for —t± < E < t± it decreases compared 
to the single layer case, d) DOS at sites B of both layers is 
equal and for —t± < E < t± it increases compared to the 
single layer case. 



III. RESULTS AND DISCUSSIONS 

A. Clean limit 

Before getting in to the discussion of the impurity 
states, we first discuss the effects of V and t± on the 
DOS in the absence of the impurity. For V = 0, and 
t± — 0, the DOS at every sites is equal and it vanishes 
linearly in E close to the Fermi energy as shown in Fig. 
[Ha). For tj_ = but y 7^ the DOS in two layers is dif- 
ferent. The overall variation of the DOS is still preserved 
but the DOS curve of upper (lower) layer shifts towards 
^ positive (negative) energy region by V/2, Fig. [D^b). So, 
the minima in the DOS for upper (lower) layer lies at 
E = \V\/2. This results in a finite density of states at 
the Fermi energy. The shift in the position of the DOS 
minima will have nontrivial influence on impurity states 
as discussed below. 

The difference in the DOS of single and bilayer 
graphene can be studied by using finite value of t±. In 
Fig- mjc-d) we show the DOS at four sites for various 
values of ij_ but fixed V = 0. The DOS at Ai,A2 is 
equal and similarly DOS at Bi, B2 is equal but those at 
A sites and B sites are not equal in the energy range 
—t± < E < t±. In this energy range, the DOS at A 
sites is smaller than that at B sites. In particular, at 
the Fermi energy, the DOS at B sites is finite but that 
at A sites is zero. These effects can be understood in 
terms of the difference in the band structure caused by 
finite t±. The band structure corresponding to sites Ai 
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FIG. 4: (Color online) For single impurity at site A\ DOS as a 
function of energy is shown. DOS at the impurity site, A\, is 
very small. DOS at the neighboring site in the same layer, B\ , 
is very high which is the signature of the impurity resonance. 
The resonance is close to the Fermi energy. Qualitatively, the 
DOS at site Ai is similar as in non impurity case. Fig. [2jc), 
but it changes at site B-2,. There is small oscillation in the 
DOS close to Fermi energy at site B2- 



FIG. 3: (Color online) LDOS in the absence of impurity at 
fixed tx = 0.1 and different V . The minimum of the DOS 
is at |V|/2. Small enhancement of the DOS occurs close to 
Fermi energy with the increase in V . We will show later that 
the position of this enhancement depends on t^. There is a 
rearrangement of DOS between A, B sites. e)and f) LDOS at 
site A\ and B\ (respectively) at fixed V = 0.5 and different 
t_L is shown. We see gap opening around the Fermi energy. 
The magnitude of the gap increases with the increase in t^. 
LDOS at site A2,B2 can be obtained by mirror inversion to 
the curves in Fig. |3je) and [3jf) about the axis oi E — 0, 
respectively. 



and A2 describe the anti-bonding states characterized by 
band gaps of ±t± which results in the decreased density 
of states compared to that of the single layer graphene 
where corresponding bands are gapless. On the other 
hand, the band structure corresponding to the sites Bi 
and B2 is gapless at the Fermi energy and has more flat 
band compared to the single layer case resulting in finite 
density of states at the Fermi level. An interesting point 
is that the DOS at B sites is also linear even though 
the dispersion relation of the electron corresponding to 
B sites can be approximated by parabolic dispersion at 
low momenta. 

Here, we also discuss the combined effect of bias and 
interlayer hoping on the DOS. In Fig. [S^a-d) we have 
shown the DOS at four sites for fixed t± = 0.1 and dif- 
ferent V. In Fig. [Slje-f) we have also shown the DOS 
at Ai and Bi for fixed V = 0.5 and different t±. Some 
new features are seen in these figures. The DOS at all 
sites gets modified compared to the unbiased case (com- 



pare with Fig. life) and (d)). There is a shift of DOS 
minimum to E = zLV/2 compared to V = case. There 
is also a gap opening up around the Fermi energy which 
increases with the increase in t±. These unique features 
of band structure will lead to non-trivial impurity states. 

B. Impurity states 

We now discuss a single impurity in the absence of 
the external bias. First we put an impurity at and 
study the DOS for various values of t±. The result of 
the actual calculation is shown in Fig. 3) The impu- 
rity has the following effects: The DOS at the impurity 
site, Ai, decreases sharply. The DOS at site Bi increases 
sharply which is the signature of a virtually bound impu- 
rity resonance state. The position of the resonance peak 
shifts slightly below the Fermi energy with the increase 
in t±. This shift of the resonance peak is expected be- 
cause the t± also affects the pole of the Greens function 
as determined by the band structure. The width of the 
resonance peak decreases with t±. This is expected be- 
cause increase in t± depresses the band DOS on site Ai, 
which suppresses the scattering rate from the generated 
impurity state. It corresponds to an increase of lifetime, 
signifying a sharp resonance. 

The effect on the DOS at A2 due to the impurity at 
Ai is very small. The overall behavior of the DOS curve 
does not change as shown in Fig. lUJc). We see a finite 
DOS at A2 at the Fermi energy which is nothing but the 
finite size effect, which is not seen in Fig. [2][c) because 
of the extended scale used to show the DOS. When we 
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UatB,,V=0.0, DOS at B, 
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(c) UatB,,V=0.O, DOSatAj (d) U at B,, V=0.0, DOS at 




FIG. 5: (Color online) For single impurity at site Bi, DOS as 
a function of energy is shown. DOS at the impurity site, Bi, is 
very small. DOS at the neighboring site in the same layer, Ai, 
is very high which is the signature of the impurity resonance. 
The resonance is close to the Fermi energy. Qualitatively, the 
DOS at site A2 is similar as in non impurity case, Fig. [2j;, but 
it changes at site B2. There is small oscillation in the DOS 
close to Fermi energy at site 52- 



zoom in and compare Fig. ^c) and Fig. we see 

small enhancement in DOS at site A2 in the presence 
of impurity at Ai (except for t± = 0). We see some 
new feature in the DOS at B2, namely the dip- hump 
structure in DOS close to the Fermi energy. We believe 
that the small change in the DOS at A2, B2 is due to the 
Friedel oscillation in the bottom layer created due to the 
impurity at the top layer. The strength of the Friedel 
oscillation in the bottom layer is very small because of 
the geometry effect. 

Now we turn to discussion on the DOS when the impu- 
rity is at site Bi, Fig. El We find that the DOS at site Bi 
is sharply reduced. At all other sites the DOS increases 
compared to the clean case, which is different from the 
result of the single impurity at Ai case. The DOS at 
Ai increases sharply which signifies the virtual bound 
state due to the impurity resonance. We notice that the 
height of the resonance peak at site Ai decreases with the 
increase in t±. Simultaneously, the resonance becomes 
more broader. This effect is just opposite to the effect 
seen in the case when the impurity is at site Ai. This 
effect can be understood in terms of the band structure 
corresponding to the site B. The increase in tj_ enhances 
the band DOS on site Bi, which increases the scattering 
rate from the generated impurity state. It corresponds 
to a decrease of lifetime, signifying a broader resonance. 
There is no much change in the overall behavior of the 
DOS at site A2 compared to the no impurity case, Fig. 
[2jc). At site B2 we see some new features, namely the 
dip-hump structure. Compared to the similar effect seen 
at i?2 due to impurity at site Ai, the structure of DOS 
in the present case is weaker. It is because the effect gets 
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FIG. 6: (Color online) LDOS for single impurity at site Ai 
with t± —0.1 fixed for different values of V. There is impurity 
resonance at site Bi. The impurity resonance is at S = V/2. 
We also see a enhancement of the DOS at site A2 assisted 
by the presence of impurity at Ai. We also see insignificant 
change in DOS at site B2- 



communicated to site B2 through Ai and A2 which will 
be weaker than the previous case where the effect was 
communicated through only one bonding length. 

Having discussed the effect in the single impurity case, 
it is much easier to understand the double impurity case. 
When there are two impurities at sites Ai and A2, the 
impurity resonance occurs at site Bi and B2. The DOS 
at Bi and B2 are identical. The height of the resonance 
peak increases with the increase in t± and the resonance 
becomes sharper and sharper. When there are two impu- 
rities at sites Bi and B2, the impurity resonance occurs 
at site Ai and ^2- The DOS at Ai and A2 are identi- 
cal. The height of the resonance peak decreases with the 
increase in t± and the peak width becomes broader. 

Now we proceed to discuss the effect of bias on the im- 
purity states in bilayer graphene. When there is a single 
impurity at site Ai , we see impurity resonance at site Bi , 
Fig. ^h) . The position of the resonance moves above the 
Fermi energy by V/2. The height of the resonance de- 
creases slightly with the increase in the V. We also see 
small satellite peak near the Fermi energy. At site A2 
also we see impurity assisted enhancement in DOS on 
the gap edge, (compare with Fig. At site B2 the 

DOS does not show much of qualitative change. 

To better understand the results, especially the satel- 
lite peak of Fig. ^h) and the position of the enhance- 
ment peak of Fig. ^c) and (d), we study the change in 
LDOS due to change in t± for fixed V = 0.5. The result 
is shown in Fig. [71 We find some new information from 
this calculation. The strength of the satellite peak, seen 
in DOS at Bi, grows with t±. We find that the position 
of the peak is determined by both V and t±. The cal- 
culation shows that the position of the satellite peak is 
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FIG. 7; (Color online) LDOS for a single impurity at site Ai, 
for V = 0.5 fixed, and with the variation in t±. 
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(c) UatB,,t^=0.1, DOSatAj (d) U at B,, t^=0.1, DOS at Bj 
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FIG. 8: (Color online) In addition to a resonance which fol- 
lows V/2, we see some interesting features at sites A2, -B2. 



This is equal to the energy gap of 



the second dip in the dispersion relation of the bilayer 
graphene in the presence of the external bias. The edge 
of the band gap seen in Fig. [3] is also dX E = Eg. Similar 
behavior is seen at site A2. The enhancement of the DOS 
seen in Fig. for V — 0.5 increases with t^. The po- 
sition of the peak also moves away from the Fermi energy 
by Eg. At site B2 new features are seen. The height of 
the small satellite peak increases with t^. The position 
of the peak also shifts symmetrically above and below 
the Fermi level by Eg. With the increase in t±^ a clear 
sign of the gap opening is seen around the Fermi energy. 

When there is a single impurity at site Bi (see Fig. [5]), 
we see different features than when there is an impurity 
at Ai . AiV — Q there is a single impurity resonance close 
to the Fermi energy. When there is an external bias, the 
resonance peak at site Ai splits in two peaks. One of 
the peak lies dX E = V/2 above the Fermi energy. Other 
remains close to the Fermi energy. With the increase 
in V, the intensity of the peak at V/2 increases and that 
close to the Fermi energy decreases. At site A2 and B2 we 
see impurity assisted enhanced DOS close to the Fermi 
energy. 

To better understand the result presented in Fig. [5] we 
perform another calculation for LDOS at fixed V — 0.5 
and different t±. The result is shown in Fig. [^1 In the 
Fig. [nja) we can see that the weight of LDOS at -E = V/2 
decreases with t±^ . This loss of weight is transferred close 
to the Fermi energy. This phenomenon can be under- 
stood in terms of the change in LDOS (at fixed V — 0.5 
in the absence of the impurity) due to the change in t^, 
(see Fig. [31[e), and[3f)), where we see that finite ij_ opens 
up a gap around the Fermi energy and the gap becomes 
more and more well defined for increased i_L. This gap 
leads to a new resonance state around the Fermi energy 
in the presence of the impurity. The shift of the weight 
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FIG. 9: (Color online) LDOS for a single impurity at site Bi, 
for V = 0.5. and with the variation in t±. 



of the LDOS close to the Fermi energy at site A2 and 
B2, can also be understood using the same logic. We 
can also see in Fig. that the LDOS at site Bi at 

E = V/2 increases with the increase in t± which explains 
the broadening effect seen in Fig. near E = V/2. 

When we have two impurities at Ai and A2 we see 
impurity resonances at Bi and B2. The position of the 
impurity resonance is again determined by V/2. Inter- 
estingly, the impurity resonance at site Bi lies above the 
Fermi energy whereas that at site B2 is below the Fermi 
energy. The DOS at Bi can be obtained by mirror re- 
flection of the DOS at B2 about the axis of E = 0. 

When we put two impurities at Bi and B2 impurity 
resonance arises at ^1 and A2. The position of the impu- 
rity resonance is again determined by V/2. The impurity 
resonance at site Ai lies above the Fermi energy whereas 
that at site A2 is below the Fermi energy. The DOS at 
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Ai can be obtained by mirror reflection of the DOS at 
A2 about the axis E — 0. 

IV. CONCLUDING REMARKS 

Before concluding, we would like to make two com- 
ments. To our knowledge two studies have been done 
which resembles to our work. Wang et. ali^^ have stud- 
ied the effects of voids in the absence of external bias. 
In particular they investigated the quantum interference 
pattern in a specific energy and no impurity induced reso- 
nance states were discussed. Similarly, the disorder prob- 
lem in biased bilayer graphene was studied by Nilsson 
and Castro Neto3 In our work we have addressed the 
issue about the sensitivity of the impurity induced reso- 
nance states to the underlying electronic band structure 
in a bilayer graphene, through the application of pressure 



or/and electric bias, which were not touched in these lit- 
eratures. 

In conclusion we study the possibility of tuning impu- 
rity states in a bilayer graphene. We systematically study 
the signature of the impurity resonance states by looking 
at the local density of states at four inequivalent sites 
of the bilayer graphene. This work should be the most 
detailed study of the impurity states in bilayer graphene. 
We have shown that in bilayer gaphene the impurity 
states can be tuned using an external bias and chang- 
ing the interlayer hopping energy. Our predictions about 
the evolution of the impurity states in bilayer graphene 
can be tested by the scanning tunnelling microscopy. 
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